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ABSTRACT 

Conventionally kinematical constraints in multibody systems are treated 
similar to geometrical constraints and are modeled by constraint reaction 
forces which are perpendicular to constraint surfaces. However, in 
reality, one may want to achieve the desired kinematical conditions by 
control forces having different directions in relation to the constraint 
surfaces. In this paper the conventional equations of motion for multibody 
systems subject to kinematical constraints are generalized by introducing 
general direction control forces. Conditions for the selections of the 
control force directions are also discussed. A redundant robotic system 
subject to prescribed end-effector motion is analyzed to illustrate the 
methods proposed. 


1. INTRODUCTION 


In many applications of multibody systems certain points are desired to 
follow prescribed paths, such as the end-effector in a robotic system. 
Such kinematical conditions are treated constraint equations to 
determine the system motion and the control forces. 

In this paper those constraints which arise from geometrical 
restrictions such as closed loops and physical guides are termed 
geometrical constraints. On the other hand, kinematical constraints are 
defined as those conditions which represent desired motions or desired 
paths of certain points or bodies. 

In the conventional methods of analysis, regardless of the 
fundamental dynamic equations (Newton-Euler, Lagrange, Kane, etc.) used, 
the constraints are modeled by constraint reaction forces which are 
perpendicular to the constraint surfaces. (See Arnold [1], Hemami and 
Weimer [2], Kamman and Huston [3], Wehage and Haug [4], Nikravesh [5], 
Kim and Vanderploeg [13], Amirouche and Jia [6].) 

However kinematical constraints do not have to be satisfied by 
constraint reaction forces, and usually have to be realized by control 
forces applied by the actuators in the system. Hence the conventional 
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solution procedure imposes an arbitrary restriction to the directions of 
the control forces. Depending on the places of the actuators in the system 
one may want to achieve the desired kinematical conditions by control 
forces having different directions in relation to the constraint surfaces. 

In this paper the conventional equations of motion are generalized by 
inroducing general direction control forces for kinematical constraints, 
that replaces the constraint force representation. And the dynamic 
equations for multibody systems subject to geometrical and kinematical 
constraints are developed. By the proposed method of solution the control 
forces and the system motion are solved simultaneously. 

This paper is divided into five sections. After the introduction, the 
second section outlines the conventional equations of motion for 
constrained multibody systems. In the third section the general direction 
control forces for kinematical constraints are introduced and the 
conditions for the control force directions are discussed. In section four 
simulations of a redundant manipulator by the conventional and the 
proposed methods are presented. Conclusions form the last section. 


2. CONVENTIONAL EQUATIONS OF NOTION 


Consider a mechanical system where qi ,... ,q n area set of generalized 
coordinates chosen for convenience to specify the configuration of the 
system. Let the system be subject to c constraints. Kane’s equations for 
an arbitrary system of particles and rigid bodies can be expressed as 
(Kamman and Huston [3], Baumgarte [7]), 


F* + F + F c = 0 


( 1 ) 


where 

p=1,...,n (2) 


F£ = Ai 


° ■ i 

*yp 


1“1,««.,C , 


F* , F and F c are the vectors of generalized inertia, external and 
constraint forces respectively. In equation (2) fi=o, 1=1,..., c are the 
constraint equations in the acceleration level, Ai are undetermined 
multipliers, and yi,...,yn are the generalized speeds of the system chosen 
for convenience as. independent linear combinations of q p . The 
tranformation between qp and y P , e.g. Euler angle derivatives and relative 
angular velocity components can be expressed as 

qh = Th P y P h,p=1 , . . . ,n ( 3 ) 

where Th P are functions of q p (Kane and Levinson [8]). 

The generalized inertia forces can be expressed in the following form 
(Huston and Passarello [9]), 
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(4) 


F* = M y + Q 

where M is the generalized mass matrix of the unconstrained system being 
functions of q P) and Q contains the quadratic velocity terms. 

The holonomic and nonholonomic constraint equations can be expressed 
in the acceleration level as below 


Bi p yp - hi i-1 , . . . ,c 

In eq. (5) B is cxn constraint matrix, and Bip and hi are in general 
functions of q P and y P . 

Note that for holonomic constraints Ch (q P ,t)=0, 


B ip = pi T hp 
and for velocity level 


nonholonomic constraints Yi (q P ,y P ,t)=0, 


Jip - 


avi 


Then, 


dy P 

using eq. 


(2), the generalized constraint forces are 


F c = B T A 


The undetermined multipliers Ai represent the restraining constraint 
forces and moments generated by the constraints at the points of 
application (Ider and Amirouche [11]). 

Equations (1) and (5) represent the governing dynamical equations. 
Combining these and making use of equations (4) and (6), we have 


M 

B T 

y ' 


F-Q ‘ 

.B 

0 

A. 


h 


The accelerations obtained from eq. (7) are then used for numerical 
integration for the time history of y P and, through the use of eq. (3), 
the generalized coordinates q P . 

Lagrange multipliers could be eliminated to reduce the equations for 
computational efficiency. To this end let C represent a nx(n-c) matrix 
which is orthogonal complement to B, obtained either by Singular value 
decomposition (Singh and Likins [12]), Zero eigenvalue method (Kamman and 
Huston [3]), Q-R decomposition (Kim and Vanderploeg [13], Amirouche and 
Jia [6]), or row equivalence transformation (Ider and Amirouche [10]). 
Premultiplying eq. (1) by C T yields 

C T (F*-F) = 0 

Combining equations (8) and (5) and utilizing eq. (4), we obtain the 
reduced equations 
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( 9 ) 


'C T M 


C T (F-Q) ' 

. B 

y = 

h 


Equations (9) and (3) form a set of 2n first order ordinary differential 
equations that can be numerically integrated to obtain the time history of 
y P and q P . 

When relative joint coordinates are selected as the generalized 
coordinates and the corresponding partial velocity vectors are developed 
using recursive multibody kinematics (Huston and Passarello [9], Ider and 
Amirouche [10]), constraint equations for joint connections are 
automatically eliminated. Hence, in this paper an open tree-like system 
represents an unconstrained system where n is the total number ofthe free 
joint degrees of freedom. 


3. CONTROL FORCES FOR KINEMATICAL CONSTRAINTS 


Now consider that a tree-like multobody system is subject to geometrical 
and kinematical constraints. Kinematical constraints represent desired 
motions or desired paths of certain points or bodies. They are the 
conditions that have to be realized by the actuators in the system. The 
desired motions could be specified at position, velocity or acceleration 
levels and could be holonomic or nonholonomic. 

Whether one uses Newton-Euler, Lagrange or Kane’s equations or other 
variations of these, in the conventional approach the constraints in the 
system are modeled by constraint reaction forces which are perpendicular 
to the constraint surfaces. In the case of geometrical constraints 
perpendicular reaction forces at the application points are generated, 
hence the above approach is necessary. On the other hand kinematical 
constraints could be achieved by a number of alternative control forces 
whose directions in the generalized space can be be selected by physical 
considerations. Therefore the conventional equations of motion should be 
generalized by considering general direction control forces for 
kinematical constraints. 

Consider c constraint equations (5), and let ci of the constraints in 
the system be geometrical and the remaining C 2 (c 2 =c-ci) be kinematical. 
The constraint matrix B and the vector of constraint force magnitudes X 
can be partitioned such that 

B = [B« t b* t ] t (io) 
and 

X = [X qT A kT ] t (11) 

where B G is a cixn matrix, B K is a C 2 matrix, A G is a ci dimensional 
vector and X K is a C 2 dimensional vector. 



Addition of control forces to the equations of motion yields 

M y + Q + B qT X g + B kT A k + A t H = F (12) 

where A is a C 2 xn control force matrix where each row represents the 
direction of the control force for each kinematical constraint in the 
generalized space, and ft Is C 2 dimensional vector of control force 
magnitudes. Now assume that the control force directions and magnitudes 
are selected such that the restraining constraint forces X G become zero. 
This leads to 

M y + Q + B qT A g + A t /l* = F (13) 

Eq. (13) can be written in the following form 


M y + Q + Z T 1/ = F 


(14) 


where Z is the augmented matrix of constraint and control force directions 
Z = [B gT A t ] t (15) 

and is the vector of constraint and control force magnitudes, 

i/ = [A gT f* T ] t (16) 

Once the control force directions A are selected by physical 
considerations eq. (14) needs to be solved together with eq. (5) to 
determine the control force magnitudes simultaneously with the generalized 
accelerations. Hence the augmented equations of motion are 



(17) 


Alternatively the equations could be reduced in a manner similar to 
Section 2. To this end, let C be a nx(n-c) matrix orthogonal complement to 
Z. Premultiplying eq. (14) by C T and augmenting with eq. (5) leads to 


C T M 

« 

C T ( F-Q ) ■ 

B 

y - 

h 


(18) 


In the case when the reduced equations are used v can be obtained from eq. 
(14) utilizing the computed accelerations, as 

V - (Z Z 1 )’ 1 Z (F-My-Q) (19) 
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Note that A should be selected such that rank Z=c, because otherwise there 
will be less than cz control forces to control C 2 kinematical conditions. 

The augmented equations have a solution if and only if the augmented 
mass matrix is nonsingular, in which case the prescribed conditions are 
achieved with the corresponding control forces. On the other hand if it is 
not physically possible to realize the kinematical constraints with the 
selected control force directions, this reveals itself as a singular ( or 
near singular) augmented mass matrix. Therefore the condition for the 
existence of solution could be expressed as follows: Directions A should 
be chosen such that a linear combination of the rows of C T M should not be 
a linear combination of the rows of B. In other words the vector space 
spanned by the rows of C T M and the vector space spanned by the rows of B 
should be non intersecting (except the zero vector). Since the dimension of 
the vector space spanned by the rows of B is c, and that of C T M is n-c, 
the possibilities for C T M are infinitely many (provided that n-c>0). Hence 
one can construct various vector spaces for C T M by different selections of 
the control force directions A. C T M that correspond to directions B is 
only one of them. 

For redundant systems, i.e. when c<n, it has been observed that one 
usually has several physically meaningful control force directions to 
satisfy the given kinematical conditions. 

In the special case when A is selected such that its rows are linear 
combinations of the rows of B K , then since 5 is the same as C in the 
conventional model, y becomes the same as the conventional case and ZV 
becomes equal to B T \. However v, may be different than Ai depending on A. 

Similarly for nonredundant systems, i.e. n-c, B is nxn, and rows of Z 
are necessarily linear combinations of the rows of B. In this case C is 
null matrix and the above procedure reduces to the conventional method 
where Vv is equal to B T \. 

It should be noted that C T M and B may form nonintersecting vector 
spaces even if C T and B do not. Hence realization of the prescribed 
motions is possible even in the extreme case when a control force 
direction is tangent to the corresponding constraint surface. This is due 
to the Inertia coupling between the generalized coordinates. 


4. SIMULATIONS OF A REDUNDANT ROBOTIC SYSTEM 


In the three link manipulator shown in Figure 1, the configuration of the 
system can be described by three generalized coordinates 0i , 02, 03. The 
generalized speeds y P are defined as 

yi = 0i , y 2 = 02 , y3 = 03 

The data used are Li =L2=L3=1 .Om, mi=30kg, m2=m3=18kg, Ii =10 kg.m 2 , 
12=13=8.64 kg.m 2 . 

The end-effector (point A) is desired to move in the horizontal 
direction with a constant velocity v A . Hence the constraint equations In 
the system are 
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Figure 1. Three link manipulator 


L i c i + I2C12 + L3C13 = v A t + 1.9088 

( 20 ) 

L 1 si + L2S12 + L3S13 = 0.9893 

where ci=cos0i, ci 2 =cos( 0 i+ 02 ), ci 3 =cos( 0 i +02+03 ) , and similarly si=sin0i, 
si2=sin(0i+02), si3=sin(0i+02+03). At the acceleration level the 

constraint equations are given by eq. (5) where B and h are 


B = 


Li Si + L2S1 2 + L3S1 3 
Li ci +L2C1 2+L3C1 3 


L2S1 2 + L3 si 3 
L2C1 2+L3C1 3 


L3 Si 3 
L 3 C 1 3 


( 21 ) 


and 


h = 


’ -Li si yi 2 -L 2 Si 2 (yi+y 2 ) 2 -Lasi3(yi+y2+y3 ) 2 ' 
. -Li ci yi 2 -L 2 Ci 2 (yi +y 2 ) 2 -L3Ci 3 (yi +y2+y3 ) 2 . 


( 22 ) 


Initially the system is at the configuration 01=60°, 02=-1O°, 
03 =-90° . The initial generalized speeds are yi=0.386 rad/sec, y2=0, 
y3=-0.9618 rad/sec, which correspond to v A = -1 m/sec. 
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First the system is simulated using the conventional method. The 
generalized constraint forces can be expressed from equations (6) and (21) 
as, 


’Ft 


Xl (L1S1+L2S12 + L3S13 ) + \2 (LlCl+L2Cl2 + L3Cl3) ' 

Fi 

= 

\i (L 2 S 12 +L 3 S 13 ) + \z (L 2 C 1 2 +L 3 C 1 3 ) 

f! 


Xl L3S13 + \Z L3C13 


In particular we wish to determine joint moments denoted as Mi, i = 1 , 2 , 3 
that would achieve the desired kinematical conditions. Ft, i = 1 ,2,3 
represent the required joint moments. It is seen from eq. (23) that all 
three joint moments are nonzero, i.e. motors are required at all three 
joints. 



Figure 2. Displacements. 

Conventional method : I.AGi, 2.A02, 3. A03 
Control force method: 4. A0i , 5.A02, 6. A03 
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Figure 3. Velocities. 

Conventional method : 1. 0i , 2. 02, 3. 03 
Control force method: 4. 0i , 5. 02, 6. 03 


The simulation is performed for 1 sec., until the end effector moves 
1m in -x direction. a0i and 0i , i=1,2,3 are plotted in Figures 2 and 3 
respectively. The joint moments Mi required to obtain the desired motion 
of the end effector as obtained by the conventional method are shown in 
Figure 4. 

Second the system is resimulated by the proposed control force 
approach. As an illustration the control force directions are selected 
such that no moments are needed at the lower joint of link 1. The 
corresponding control force directions are 


A 



1 0 ' 

0 1 . 


(24) 


Note that since there are no geometrical constraints in this system, the 
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Time (see) 



Figure 4. Joint moments. 

Conventional method : 1. Hi, 2. M 2 , 3. M 3 
Control force method: 4. M 2 , 5. M 3 


matrix Z is identical to A, and the vector v is identical to /u . The 
control forces Z T v are 



0 


0 


0 

Z T V = Vi 

1 

+ V 2 

0 

= 

Vi 


. 0 . 


.1 . 


.v 2 


(25) 


Hence, in this case, the required joint moments for the prescribed end 
effector motion are Mi=0, M 2 = vi, M3=V2. 

The augmented mass matrix was observed to be full rank as expected. 
9i and 0i , i = 1 ,2,3 for a simulation of 1 sec. are plotted in Figures 2 and 
3. The joint moments are plotted in Figure 4. Note that in this case a 
motor is not needed at the lower joint of link 1. 
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2 


1 



Initial configuration 


Figure 5. Final configurations: 1. Conventional 
method; 2. Control force method. 


The system’s motion differ in the above two approaches although in 
both cases the end effector performs the same motion. The configurations 
at t=1 sec. corresponding to the conventional model and the control force 
model are shown in Figure 5. 


5. CONCLUSIONS 


This paper presented a general procedure for the dynamic modeling of 
multibody systems subject to kinematical constraints. General direction 
control forces have been introduced that replace the conventional 
constraint reaction forces, hence increasing the ways of realization of 
the prescribed motions. It is shown that the possible conrol force 
directions are more than one, and the criteria for the existence of 
solution have been presented. 

The method proposed in this paper involves selecting the control 
force directions in the generalized space by physical considerations, and 
then solving their magnitudes simultaneously with the corresponding system 
motion. As a result one can design alternative control forces that can be 
applied by the actuators in the system. The method is expected to be 
especially useful to control the extra degrees of freedom in systems that 
have joint flexibility or joint clearance. 
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